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I.    INTRODUCTION 
Jet  Vane  TVC  Systems 

Thrust  vector  control  (TVC)  systems  offer  means  of 
flight  vehicle  trajectory  control  that  are  virtually  inde- 
pendent of  external  forces.   Such  a  capability  is  frequently 
required  for  tactical  missiles,  as  well  as  spacecraft  launch 
vehicles,  when  the  relative  flow  past  external  lifting  sur- 
faces is  insufficient  to  generate  the  necessary  control 
forces.   This  commonly  occurs  during  low-speed  operations, 
such  as  at  launch  or  during  hovering  flight.   High  angle  of 
attack  flight  may  also  lead  to  situations  in  which  conven- 
tional lifting  surfaces  are  inadequate.   In  addition,  there 
are  occasions  when  external  steering  devices  are  infeasible 
from  a  design  point  of  view,  such  as  for  tube-launched 
devices. 

Several  methods  of  TVC  have  been  developed  and  applied 
to  operational  and  experimental  vehicles.   These  include 
movable  nozzles,  internal  fluid  injection  (secondary  in- 
jection), and  mechanical  jet  deflection  systems.   Jet  vane 
systems  fall  in  the  latter  category  and  they  tend  to  be 
favored  for  volume-limited  applications  requiring  relatively 
low  actuation  torques,  large  thrust  deflection  angles,  and 
rapid  response.   Jet  vanes  may  also  be  used  with  relative 


ease  to  generate  roll  torques.   The  application  of  Jet  vane 
TVC  dates  back  to  the  rockets  designed  by  Goddard,  and  has 
extended  to  the  Redstone,  Sergeant,  Talos,  Pershing,  and 
Algo  II  and  III  motors  CI, 2],  as  well  as  several  instal- 
lations in  smaller  tactical  rockets. 

Of  course,  there  are  disadvantages  accompanying  the 
selection  of  jet  vanes  for  TVC  purposes.   These  include 
thrust  losses  on  the  order  of  3-5%  with  undeflected  vanes 
C2].   In  addition,  the  attainment  of  relatively  high  thrust 
deflection  angles  may  lead  to  axial  thrust  losses  of  the 
same  order  of  magnitude  as  the  resulting  side  force.   How- 
ever, the  chief  problem  associated  with  the  use  of  jet  vanes 
is  the  large  thermal  loading  that  they  experience  as  they 
are  required  to  operate  in  hot,  high-speed,  particle-laden 
flows.   This  problem  leads  to  design  limitations  so  that  jet 
vanes  are  often  restricted  to  short-duration  use  in  motors 
with  low-temperature  non-metalized  propellants. 

The  aerodynamic  (side-force  producing)  characteristics 
of  jet  vanes  may  be  calculated  with  fair  certainty  on  the 
basis  of  inviscid  flow  theory  with  suitable  corrections  for 
viscous  effects  C2, 31.   On  the  other  hand,  difficulties  that 
stem  from  the  severity  of  the  Jet-vane  thermal  environment 
have  led  to  design  practices  that  are  based  largely  upon 
past  experience  and  cut-and-try  methods.   Over-design  is 
therefore  inevitable,  with  virtually  no  capacity  for  design 
optimization.   In  order  to  exploit  the  several  advantages  of 
Jet  vane  TVC   systems,  therefore,  it  has  become  necessary  to 


build  reliable  data  bases  and,  to  the  extent  possible,  at- 
tain a  fundamental  understanding  of  the  heat  transfer  char- 
acteristics of  such  systems. 

This  need  has  led  to  the  work  undertaken  at  the  Naval 
Postgraduate  School  ( NPS )  in  support  of  a  larger  program  at 
the  Naval  Weapons  Center  (NWC).   Previous  investigations  at 
NPS  have  included  applications  of  computational  fluid  dyna- 
mics (CFD)  [4,5]  and  wind-tunnel  tests  using  infra-red  ther- 
mography [6D.   A  summary  of  the  results  of  these  studies, 
together  with  an  overview  of  relevant  previous  works,  are 
contained  in  Ref.  C7],   Work  in  the  area  of  CFD  has  con- 
tinued C8]  and  recent  results  [9]  have  shown  that  dynamic 
simulation  methods  hold  promise  in  further  identifying  the 
dominant  factors  affecting  the  thermal  characteristics  of 
jet  vanes.   This  latter  area  of  study  has  advanced  to  the 
application  of  parametric  system  identification,  and  it  is 
the  results  of  these  efforts  that  are  the  main  subject  of 
this  report. 
The  Jet  Vane  Thermofluid  Environment 

In  the  design  of  a  jet  vane  system,  the  integrity  of 
the  vanes  themselves  must  be  guaranteed  over  the  specified 
work  cycle.   Although  this  is  a  serious  challenge,  designers 
must  also  consider  the  behavior  of  the  vanes  and  supporting 
structure  during  transient  events.   Upon  motor  ignition, 
localized  jet  vane  temperatures  may  rise  to  near-stagnation 
values  within  a  few  seconds.   Temperatures  in  the  vane  at- 
tachment device  will  rapidly  follow  this  rise,  and  severe 


thermal  stresses  may  develop  due  to  the  proximity  of  a  rela- 
tively cool  supporting  structure.   These  and  other  design 
aspects  can  only  be  addressed  with  precision  if  there  is  a 
good  understanding  of  the  convective  heat  transfer  process 
that  gives  rise  to  the  energy  transfer  from  the  flowing 
gases  to  the  the  vane.   Put  another  way,  the  application  of 
complex  computer  codes  for  thermal  conduction  in  the  vane 
and  supporting  structure  can  only  follow  the  specification 
of  the  convective  boundary  conditions. 

With  respect  to  these  conditions,  the  problem  is  even 
more  complex.   The  vane  is  immersed  in  a  flow  field  that  is, 
if  generally  described,  compressible,  turbulent,  multi-com- 
ponent (and  possibly  multi-phase),  three-dimensional,  and 
unsteady,  with  variable  properties  and  nonlinear  and  time- 
variant  boundary  conditions.   Even  if  taken  one  at  a  time, 
these  complexities  present  problems  that  are  beyond  the 
state-of-the-art  for  exact  solution.   In  addition,  the  over- 
all flow  field  will  contain  intersecting  and  impinging  shock 
waves  that  give  rise  to  discontinuous  events  and  further 
complication  of  the  boundary  conditions.   The  presence  of 
various  protuber-ances  only  serves  to  exacerbate  these 
difficulties. 

To  further  define  the  problem,  it  may  be  of  use  to 
consider  the  levels  of  heat  transfer  that  might  be  expected 
from  the  point  of  view  of  "simple"  convection  from  a  super- 
sonic flow  to  a  cooled  wall.   The  convective  heat  transfer 
coefficient  (h)  may  be  described  in  terms  of  the  nondimen- 


sional  Stanton  number  (St)  as  follows: 

h  =  (pVc  )  St 

where   p,  V,  and  c   are  the  density,  velocity,  and  specific 

P 

heat  at  constant  pressure  of  the  flowing  gas.   On  the  as- 
sumption that  the  gas  behaves  ideally,  this  expression  may 
be  written 

k   P 
h  =  St  <-*-  z-  g  I   ) 
k-1  T     sp 

where  P  and  T  are  local  gas  pressure  and  temperature  and  I 

sp 

is  the  specific  impulse.   Consider,  for  instance,  the  condi- 
tions at  the  exit  of  a  rocket  nozzle  flowing  at  Mach  3. 
With  a  ratio  of  specific  heats  assumed  to  be  k  =  1.2,  we 
have,  approximately: 

P 

h  =  St  (0.25  ^9  I   > 
T      sp 
o 

where  the  subscript  (  )  refers  to  stagnation  conditions.   If 

o 

a  typical  case  is  taken  to  be  given  by  P   =  100  bar,  T   = 

o  o 

2560  K,  and  I    =  250  s,  then,  in  round  numbers,  the  heat 
sp 

transfer  coefficient,  in  units  of  (kW/maK),  is  h  =  2000  St. 

In  rocket  engine  nozzles,  a  typical  value  for  the 
Stanton  number  is  about  0.002  C10].   Thus,  with  the  result 
above,  the  expected  value  of  the  heat  transfer  coefficient 
is  on  the  order  of  4  kW/maK.   Such  a  value  should  not  be 
taken  as  conservative,  since  actual  motor  conditions  may  be 
more  severe  and  local  regions  (e.g.  stagnation  points  or 
regions  of  flow  reversal)  may  experience  much  higher  values. 
By  way  of  comparison,  values  of  h  such  as  these  are  much 


larger  than  those  usually  considered  for  convective  heat 
transfer  and,  in  fact,  they  are  more  typical  of  those  real- 
ized in  phase-change  processes  [11]. 

II.    MODELING  AND  SIMULATION 
Background 

A  general  goal  of  jet  vane  heat  transfer  studies  is  to 
develop  a  capability  by  which  the  transient  thermal  behavior 
of  the  vane  may  be  predicted  with  confidence.   An  under- 
standing of  the  energy  transport  processes  occurring  at  the 
boundaries  of  the  vane  is,  of  course,  essential  to  the 
achievement  of  this  goal  but,  because  of  the  complexity  of 
the  thermofluid  environment  in  which  the  vane  operates,  CFD 
methods  cannot  be  expected  to  yield  a  complete  and  compre- 
hensive predictive  capability.   In  addition,  the  results  of 
CFD  modeling  must  invevitably  be  simplified  if  they  are  to 
be  adapted  to  vane  design  purposes  and  the  prediction  of 
system  performance.   In  order  to  complement  the  CFD  studies, 
therefore,  a  program  has  been  initiated  to  approach  the 
problem  from  a  less-general  but  more-practical  point  of 
view.   Underlying  this  approach  is  an  acceptance  of  the  fact 
that  details  of  the  energy  transport  processes  will  never  be 
fully  defined.   On  the  other  hand,  there  is  also  the  recog- 
nition that  the  cumulative  or  "lumped"  effect  of  all  the 
various  complications  mentioned  above  is  to  transfer  energy 
to  and  from  the  vane.   These  boundary  processes  lead  to 
further  energy  flow  within  the  vane  that  ultimately  gives 


rise  to  a  vane  temperature  distribution  that  varies  with 
time. 

This  transient  behavior  should  reflect  the  nature  of 
the  boundary  heat  transfer  processes  that  drive  it,  as  well 
as  the  thermal  impact  of  various  vane  design  parameters 
(thermal  conductivity,  heat  capacity,  density,  geometry, 
etc. ).   The  transient  response  is  the  "signature"  of  the 
combined  effects  of  the  thermal  environment  and  the  vane 
configuration.   This  premise  is  the  basis  for  justification 
of  the  modeling  and  simulation  (M&S)  study  described  in 
subsequent  sections  of  this  report. 

The  goal  of  the  M&S  study  has  been  to  develop  a  dynamic 
model  of  vane  heating  and  cooling.   If  such  a  model  can  be 
constructed,  even  though  it  is  approximate,  then  significant 
benefits  will  accrue.   For  instance,  measured  vane  tempera- 
ture histories  can  be  used  for  deductive  purposes  such  as 
parametric  system  identification  (PSI)  whereby  the  model  is 
used  to  deduce  what  would  have  had  to  be  true  in  order  for 
the  observed  temperature  history  to  have  occurred.   Having 
estimated  the  values  of  appropriate  parameters,  local  tem- 
peratures can  be  deduced  at  points  in  the  vane  structure 
where  conditions  (accessibility,  sensor  survival,  etc. ) 
preclude  measurement.   In  addition,  a  trustworthy  predictive 
model  would:  <1)  reduce  the  occasions  in  which  testing  is 
required  to  verify  the  effects  of  design  changes,  (2)  permit 
the  extension  of  sub-scale  test  data  to  predict  full-scale 
jet  vane  performance,  and  (3)  indicate  the  design  directions 
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most  likely  to  lead  to  optimum  system  performance. 
Approach 

The  M&S  approach  may  be  thought  of  as  a  vast  simplifi- 
cation of  a  numerical  model  that,  were  it  possible  to  con- 
struct, would  lead  to  a  successful  predictive  capability. 
Whereas  such  a  comprehensive  model  would  treat  the  flow 
environment  of  the  vane  and  the  vane  itself  in  fine  numeri- 
cal detail,  the  M&S  method  assumes  that  sufficient  accuracy 
may  be  obtained  if  the  flow  and  the  vane  are  made  up  of  a 
relatively  few  thermal  parts.   An  Important  product  of  the 
M&S  approach  is  the  prediction  of  the  time  history  of  the 
vane  temperature  --  in  fact,  this  is  what  provides  the  vital 
clues  that  are  used  to  estimate  the  values  of  the  model 
parameters. 

The  system  identification  procedures  used  in  this  study 
are  those  familiar  to  the  discipline  of  automatic  controls. 
If  the  governing  mathematical  relationships  can  be  cast  in  a 
transfer-function  or  other  systematic  form  (not  necessarily 
linear),  then  experimental  observations  can  be  related  to 
the  various  parameters  of  the  model  by  system  identification 
methods.   In  order  for  this  procedure  to  be  of  value,  how- 
ever, the  afore-mentioned  parameters  must  be  related  to  the 
physical  quantities  affecting  the  thermal  transport  process. 
The  establishment  of  these  relationships  is  described  in  the 
next  section  of  this  report. 


III.    MODEL  DISCRETIZATION 
The  model  as  presently  configured  considers  only  two 
thermal  energy  input  processes:  forced  convection  at  the 
vane  surface  and  stagnation-point  heat  transfer  at  the 
leading  edge.  These  processes  are  driven  by  the  flow  stag- 
nation temperature  which  is,  in  turn,  derived  from  the 
measured  thrust  levels  by  means  of  the  rocket-motor  bal- 
listic characteristics.  The  rise  and  fall  of  the  stagnation 
temperature  is  assumed  to  follow  that  of  the  thrust  without 
significant  dynamic  offset.    In  other  words,  the  generic 
product  of  the  model  is  a  series  of  functions  that  give 
local  vane  temperatures  in  response  to  thrust.    (If  these 
functions  were  linear,  they  would  have  direct  transfer- 
function  counterparts. ) 

Several  vane  discretization  schemes  have  been  tested, 
and  it  has  been  found  that  four  "lumps"  are  sufficient  to 
indicate  the  general  thermal  behavior  of  the  vane:  vane  tip, 
vane  body,  shaft,  and  mount.   This  configuration  leads  to  a 
simplicity  that  is  a  necessary  feature  of  the  M&S  approach, 
and  results  to  be  described  later  justify  the  use  of  such 
simple  models  as  design  tools. 
Estimation  of  Geometrical  and  Thermal  Properties 

As  a  preliminary  step  in  the  study,  the  actual  vane 
design  of  interest  to  the  Naval  Weapons  Center  was  confi- 
gured as  a  collection  of  lumps  with  geometries  suitable  for 
estimation  of  thermal  conduction  properties.   Figure  1  il- 
lustrates the  vane  geometry  that  has  been  modeled  (the  vane 
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Figure  1.   Naval  Weapons  Center  jet  vane  configuration. 
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is  constructed  of  lOV.Cu/W,  and  the  mount  material  is  steel). 

From  this  design,  the  fictitious  thermal  vane  has  been 
hypothesized  as  consisting  of  three  rectangular  solids:  tip, 
fin,  and  shaft.   This  discretized  vane  is  shown  in  Fig.  2. 
The  tip  of  the  vane  has  been  separately  identified  in  order 
to  account  for  the  stagnation  properties  of  the  thermal  con- 
vection near  the  vane  leading  edge.    This  portion  has  been 
arbitrarily  sized  so  as  to  have  a  chord  length  of  10*/.  of  the 
total  vane  chord  of  3. 75  inches.   With  this  length  fixed, 
the  discretized  tip  is  shaped  so  as  to  have  the  same  lateral 
area  (chord  x  thickness)  as  that  of  the  actual  tapered  and 
rounded  tip. 

The  remainder  of  the  discretized  fin  is  assumed  to  be 
subjected  to  thermal  convection  of  a  turbulent  boundary 
layer  type.   It  is  lumped  into  a  rectangular  solid  of  thick- 
ness equal  to  the  average  value  of  the  tapered  fin  (with  tip 
removed)  and  span  equal  to  the  mean  span  of  the  actual  fin. 
The  remaining  length  dimension  is  set  by  equating  the  volume 
of  the  rectangular  solid  to  that  of  the  fin  (less  the  tip 
portion).   The  vane  shaft  is  similarly  "molded"  into  a  rec- 
tangular solid,  and  this  element  is  assumed  to  be  subjected 
only  to  conduction  (and  radiation,  were  it  included)  heat 
transfer. 

The  vane  mount  (not  shown  in  Figs.  1  &  2)  is  a  rela- 
tively massive  and  complex  structure.   In  order  to  pursue 
the  basic  feasibility  of  the  method,  no  attempt  was  made  to 
model  the  thermal  resistance  of  this  component  in  detail. 
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Figure  2.   Discretized  version  of  NWC  jet  vane. 
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Instead,  representative  sets  of  values  (conduction  length 
and  cross-section)  for  the  main  mount  components  were  esti- 
mated using  a  rough  scaling  procedure.   The  appropriate 
thermal  resistance  of  the  mount  was  then  computed  from  an 
analysis  of  the  analogous  electrical  circuit.    In  prelimi- 
nary studies  the  mount  thermal  resistance  was  considered  to 
be  an  "adjustable  parameter"  for  use  in  seeking  agreement 
with  the  available  experimental  evidence.    (Determination  of 
the  correct  value  of  this  quantity  is  one  of  the  goals  of 
the  system  identification  procedure  described  later. ) 

Each  of  the  components  described  above  was  assigned  a 
thermal  node  located  at  it's  mass-center  which,  in  turn,  is 
the  assumed  location  of  the  energy  storage  associated  with 
the  entire  mass.   The  capacity  to  store  energy  was  calcu- 
lated in  the  usual  way  for  each  lump.   For  the  mount,  the 
thermal  capa-city  was  assumed  to  be  infinite  --  thus  the 
temperature  of  this  node  remained  constant  at  ambient 
(ground)  temperature  during  the  simulation. 

It  should  be  emphasized  that  at  this  point  in  the  study 
the  discretization  rationale  is  quite  arbitrary  --  the  goal 
has  been  to  obtain  adequate  agreement  with  test  results 
using  the  minimum  number  of  thermal  components.   Refinement 
of  the  model  can  easily  include  the  division  of  the  vane 
structure  into  a  larger  number  of  smaller  "lumps. "   Although 
this  would  permit  the  estimation  of  temperatures  at  more 
nodes  and  interfaces,  the  accuracy  of  these  estimates  would 
be  no  greater  than  that  associated  with  the  "minimum  node" 
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model  described  here.   Such  added  complication  of  the  vane 

nodal  distribution  would  only  be  Justifiable  in  conjunction 

with  a  more-detailed  description  of  the  thermal  boundary 

conditions,  and  this  would  add  even  more  uncertainty  to  the 

model.    Results  described  later  indicate  that  such 

complication  is  not  necessary. 

Development  of  the  Governing  Equations 

With  the  thermal  network  and  associated  properties  thus 

defined,  the  governing  equations  were  formulated  by  means  of 

an  energy  balance  at  each  node.    In  general,  this  balance 

reads  "rate  of  heat  flow  in  =  rate  of  heat  flow  out  *■    rate 

of  energy  storage. "   For  the  fin  node  (node  2  in  Fig. 3),  for 

instance : 

T__  -  T_    TJ  -  T^    T_  -  T_    _ 
R2 2    _i 2  _  _2 3    C  sT 

R^^  *»  .  ^-»         R^« 

F2  i2  23 

Here  the  symbol  R  is  used  to  denote  thermal  resistance  -- 
1/hA  for  convection  to  the  fin  and  L/kA  for  conduction 
between  the  nodes  --  and  C  is  thermal  capacitance  (mass  x 
heat  capacity).   The  letter  s  is,  as  usual,  the  symbol  for 
the  Laplace  variable.   Calculation  of  the  recovery  temper- 
ature T^„  is  discussed  in  the  next  section. 
R2 

The  temperature  at  the  node  may  be  expressed  explicitly 
in  terms  of  the  surrounding  nodal  temperatures  as  follows: 

2     1  .  ,2B   Rp2    R±2    R23 

where  the  nodal  resistance  R  _   is  given  by: 

n2 


Rn2    RF2    Ri2    R23 
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NOMINAL  VALUES 


RESISTANCES 

CAPACITANCES 

(K/W) 

(J/K) 

Fl    «  20.3 

NODE  1     *  0.11 

F2    *      2.3 

NOOE  2    *  2.2 

1  i    »      0.35 

NODE  3    *   1 . 0 

i  2  «      1.6 

23  «      2.6 

3G  «    10.3 

Figure    3.      Nodal   configuration  and   table   of   estimated 
values. 
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and  the  nodal  time  constant  is  given  by 

t-  =  R  _  C 
2     n2   2 

Nodes  of  the  model  that  do  not  provide  an  energy  storage 
function  are  also  absent  the  associated  time  constant.   Thus 
the  temperature  at  the  internal  node  shown  at  the  interface 
between  the  tip  and  fin  of  the  present  model  (node  i  in  Fig. 
3)  is  expressed  in  terms  of  the  temperature  of  the  surround- 
ing nodes  with  only  thermal  resistances  as  parameters. 

A  heat  balance  at  each  of  the  remaining  nodes  yields 
relationships  similar  to  those  given  above,  so  sufficient 
information  is  available  to  determine  each  nodal  temperature 
as  a  function  of  time.    Initially,  the  system  is  assumed  to 
be  at  thermal  eguilibrium  with  the  environment  so  that  the 
nodal  temperatures  are  all  equal  to  the  ambient  value.   With 
the  firing  of  the  rocket  motor  this  equilibrium  is  disturbed 
and  the  two  recovery  temperatures  provide  inputs  to  the 
ensuing  heat  transfer  process.    In  preliminary  studies  the 
Digital  Simulation  Language  (DSL)  has  been  used  to  automati- 
cally accomplish  the  repeated  integration  steps  required  to 
determine  the  thermal  response  of  those  nodes  possessing 
storage  capacities. 

The  calculation  method,  together  with  preliminary 
results,  are  described  in  later  sections.   Before  this, 
however,  it  is  necessary  to  estimate  the  driving  functions 
of  the  entire  process  --  the  thermal  inputs  at  the  tip  and 
fin  surfaces.   This  is  the  subject  of  the  next  chapter. 
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IV.    ESTIMATES  OF  SURFACE  HEAT  TRANSFER  RATES 
As  described  in  the  previous  section,  the  forcing  func- 
tions of  particular  interest  are  those  due  to  heating  at  the 
vane  tip  and  along  the  downstream  surface  of  the  vane.    (The 
cooling  effects  of  radiation  and  ablation,  both  clearly  im- 
portant in  the  actual  situation,  are  omitted  from  discussion 
for  the  time  being. )   Most  of  the  heat  transfer  computations 
outlined  below  are  derived  from  the  analytical  methods  sug- 
gested in  the  works  edited  by  C.C.     Lin  C12]  and  in  the  AGARD 
Monograph  authored  by  Ziebland  and  Parkinson  CIO]. 

A  particular  feature  of  these  high  speed  flows  is  the 
large  differences  in  temperature  that  the  gas  experiences  in 
decelerating  near  the  body  surface.    In  such  cases  it  is 
necessary  to  account  for  the  temperature  dependency  of  gas 
properties.   A  most-useful  simplification  for  this  purpose 
is  that  proposed  by  Diessler  and  others  C 12,  p. 304]  in  which 
the  Prandtl  number  and  specific  heat  are  considered  to  be 
constant  inasmuch  as  their  variations  with  temperature  are 
of  a  lower  order  of  magnitude  than  those  of  the  other  gas 
properties  (viscosity  and  thermal  conductivity).   Thus,  with 
this  assumption,  the  quantity  Pr/c   =  constant  and  a  separ- 
ate estimate  of  the  gas  viscosity  leads  to  the  thermal 
conductivity  for  a  given  Prandtl  number  and  specific  heat. 

For  the  rough  calculations  used  in  this  study,  the 
Prandtl  number  has  been  estimated  using  the  Eucken  formula 
C2,  p. 1393,  Pr  =  4k/(9k-5),  where  k  is  the  ratio  of  specific 

heats.   With  the  gas  constant,  R  ,  the  constant-pressure 

g 


la 


specific  heat  is  given  by  c   =  R  k/(k-l). 

p     g 

In  this  work,  the  standard  Sutherland- type  formula  has 
been  adopted  for  determining  viscosity: 

ref      ref 

Following  the  recommendation  of  [2,  p. 93,  a  value  of  n  =  0. 7 

has  been  used.    In  addition,  the  reference  viscosity  has 

been  taken  to  be  u    _  =  4.0xl0~   N-s/m   at  T    =  1000  K. 

ref  ref 

Stagnation  Point  Heat  Transfer 

In  the  analysis  of  stagnation  point  heating,  the  solu- 
tion to  the  boundary  layer  equation  requires  an  estimate  of 
the  local  fluid  acceleration  in  the  vicinity  of  that  point. 
In  the  case  of  supersonic  flow,  this  can  be  approximated 
upon  the  assumption  that  Newtonian  flow  prevails  between  the 
bow  wave 

bow  wave  and  the  body  and,   (see  van  Driest  CI,  p. 366]), 
fl  =  (du  /dx)    =  (U/D)  Ca<P  /P   )(T   /T  )  ]°'5 

OD         X=0  09     Oy       OflD     00 

Here  the  subscript  ( a» )  refers  to  freestream  conditions,   (o) 
denotes  stagnation  conditions,  and  (y)  is  for  conditions 
downstream  of  the  normal  shock.   The  pressure  and  tempera- 
ture ratios  in  the  above  expressions  are  known  functions  of 
the  freestream  Mach  number. 

With  fi  thus  defined,  the  Stanton  number  may  be  computed 
from  the  following: 

0.  5   -0.  6    -0.  5 
St   =  0.57  (flD/U)    Pr      Re 
o  o       o 
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Turbulent  Boundary  Layer  Convection 

In  high  speed  compressible  flows,  the  analysis  of  ther- 
mal processes  is  complicated  by  the  fact  that  considerable 
compression  and  viscous  dissipation  follow  from  the  deceler- 
ations occuring  in  the  boundary  layer.   This  leads  to  tem- 
peratures within  the  boundary  layer  that  are  in  excess  of 
that  of  the  freestream,  and  the  driving  temperature  for  heat 

transfer  is  the  so-called  recovery  temperature,  T  .   Thus  a 

R 

recovery  factor  is  defined  as  follows: 

Tn  -  T 
R ® 

r    T   -  T 

O       co 

Fortunately,  it  has  been  found  by  many  investigators  that 
the  recovery  factor  may  be  related  to  the  Prandtl  number  in 
a  simple  way  that  is  adequate  for  most  purposes.   For  a  tur- 
bulent boundary  layer  (assumed  here),  the  relationship  is 
r  =  Pr 

A  more-difficult  problem  arises  from  the  dependency  of 
gas  properties  upon  temperature  since  in  most  cases  the 
appropriate  reference  temperature  for  this  calculation 
depends  upon  the  wall  temperature.   This,  in  turn,  depends 
upon  the  reference  temperature.   The  reference  temperature 
is  typically  defined  as  [2] 

T   _  =  0.  5  T   ,  ,  t  0.  28  T   +  0.  22  T^ 
ref         wall  co  R 

The  present  calculations  are  based  upon  the  assumption  that 
the  wall  temperature  is  adequately  represented  for  this  pur- 
pose by  the  mean  of  the  recovery  and  ambient  temperatures. 
(An  evaluation  of  the  validity  of  this  assumption  requires 
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an  iterative  time-dependent  calculation.    In  the  system 
identification  work  to  be  described,  however,  the  convective 
resistance  is  treated  as  a  parameter  to  be  identified.   Thus 
the  experimental  data  lead  to  the  deduction  of  the  "effec- 
tive" film  coefficient  and  the  detailed  analyses  described 
in  this  section  may  not  be  necessary  for  design  purposes. ) 

Given  that  the  reference  temperature  for  gas  properties 
is  adequately  described,  the  gas  viscosity  and  thermal  con- 
ductivity may  be  estimated  as  above.   The  Stanton  number  for 

turbulent  compressible  flow  is  then  given  by  C23: 

~  ~^^r     n    -0. 67  _  -0. 2 
St  =  0. 0296  Pr        Re 

From  the  expressions  given  above  for  the  Stanton  num- 
bers, Nusselt  numbers  and  thermal  resistances  can  be  calcu- 
lated for  the  model  components  affected  by  stagnation  and 
boundary  layer  heat  transfer  processes.    It  should  be  noted 
that  the  thermal  resistances  are  scale-dependent:  stagnation 

point  thermal  resistance  decreases  as  (scale)   ",  and  turbu- 
lent boundary  layer  thermal  resistance  is  propportional  to 

-1.2 
(scale)      .   The  simulations  conducted  in  this  study  have 

been  of  a  1/4-scale  vane  in  order  to  provide  a  comparison 

with  available  NWC  test  data. 

Input  Modeling 

The  information  necessary  to  provide  values  for  the 

thermal  resistance  includes  the  freestream  characteristics 

at  the  location  of  the  jet  vane  --  Mach  number,  stagnation 

temperature,  and  stagnation  pressure.   The  concept  of  the 

present  model  is  that  these  quantities  are  determined  from 
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the  ballistic  characteristics  of  the  rocket  motor.   Accord- 
ingly, inputs  required  for  the  simulation  include  the  motor 
chamber  pressure,  thrust,  and  characteristic  velocity,  and 
the  discharge  coefficient  and  pressure  ratio  of  the  nozzle. 
In  addition,  necessary  propellant  gas  properties  include  the 
ratio  of  specific  heats  and  the  gas  molecular  weight  or  gas 
constant.   Experimental  values  are  preferred,  to  the  extent 
that  they  are  available,  but  a  considerable  body  of  theory 
exists  if  analytical  estimates  are  necessary. 

From  these  quantities,  the  stagnation  temperature  may 
be  calculated  as  follows: 

T   =  ICJ  r  c*]2/R 
o      d  g 

where: 

C.   =  nozzle  discharge  coefficient  (.934) 
d 

r    =  Ck(2/k.l)(k+1>/<k-1)]°-5  (0.65) 

* 
c    =  characteristic  velocity  (1512  m/s) 

2   2 
R    =  gas  constant  (318.5  m  /s  -K) 

g 

k    =  ratio  of  specific  heats  (1.21) 
The  numbers  given  in  parentheses  are  those  presently  in  use, 
and  lead  to  a  stagnation  temperature  of  2650  K. 

From  the  motor  nozzle  pressure  ratio  Mach  number  at  the 

exit  (assumed  to  be  that  at  the  vane)  is  computed  from 

2  ( k-1 ) /k 

M         =     C2/(k-l)]     [<P       /P    >  -    1] 

OD  Q<D  CD 

and  selecting  a  nozzle  pressure  ratio  of  P   /P   =  186,  this 

OCD     OD 

expression  gives  M   =3. 75.    It  is  important  to  note  that 

OD 

the  motor  chamber  pressure  must  also  be  provided  because  the 
heat  transfer  calculations  for  the  vane  require  the  density 


22 


and  hence  the  pressure  in  the  freestream  at  the  vane  loca- 
tion --  in  the  calculations  presented  here,  the  value  of 
P    =  15. 76  MPa  has  been  used. 

Qflj 

For  ease  of  reference,  the  following  table  of  values 

describes  the  heat  transfer  quantities  obtained  from  the 

inputs  and  computed  quantities  described  above. 

Stagnation     Turbulent  Boundary 

Point  Layer 


Stanton  number  8.53x10 
Nusselt  number  43. 1 

Film  coefficient,  W/m2-K  6. 46x10' 
Thermal  resistance,  K/W         5.69 


3.  01x10 

438 

2. 25x10 

0.  634 
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V.  SIMULATION  RESULTS 
The  model  and  governing  equations,  described  above, 
have  been  coded  in  the  Digital  Simulation  Language  (DSL). 
This  IBM  software  product  is  a  FORTRAN-based  language  with 
the  particular  advantage  that  a  wide  variety  of  integration 
subroutines  are  imbedded  and  directly  available  to  the  user. 
Such  integrations  are  necessary  at  each  of  the  energy  stor- 
age nodes,  and  for  this  purpose  a  5-th  order  Runge-Kutta 
method  incorporating  a  variable  step  size  was  selected.   In 
addition  to  the  vane  conduction  model  described  above,  the 
code  incorporates  energy  flow  calculations  in  order  to  veri- 
fy that  the  model  does  indeed  account  for  the  disposition  of 
all  heat  transferred  to  or  from  the  vane,  as  well  as  that 
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stored  at  the  various  nodes. 

The  simulation  code  has  been  written  so  that  the  dri- 
ving input  is  the  thrust  of  the  rocket  motor,  as  provided  by 
test  data.   In  the  case  reported  here,  the  point-by-point 
thrust  vs.  time  measurements  are  approximated  by  a  ramp  in- 
put from  0  to  2335  N  (525  lb)  in  a  period  of  0.5s,  held  at 
this  level  for  2. 5s,  and  ramped  back  to  zero  thrust  in  a 
further  period  of  0. 5s.   Although  the  actual  thrust  data 
points  could  be  used,  the  ramp-up/ramp-down  closely  approxi- 
mates the  thrust  schedule  and  provides  a  useful  input  for 
evaluating  the  sensitivity  of  the  thermal  response  of  the 
vane  to  various  parameters  of  the  model. 

An  important  aspect  of  the  physical  event  of  rocket 
firing  is  that  the  convective  resistances  to  vane  heat 
transfer  (at  the  stagnation  point  and  in  the  boundary  layer) 
are  coupled  to  the  presence  or  absence  of  flow  past  the 
vane.   Thus  these  resistances  are  initially  very  large, 
decrease  rapidly  to  plateau  values  as  the  motor  reaches  full 
thrust,  and  increase  again  during  burnout.   This  behavior 
has  been  modeled  by  postulating  that  the  film  coefficients 
begin  and  end  at  IX  of  the  full-flow  values  and  follow  the 
ramp-up,  plateau,  ramp-down  profile  of  the  thrust  schedule. 

For  comparison  with  the  simulation  results,  the  vane 
temperature  data  were  likewise  converted  to  a  continuous 
record  by  means  of  the  transfer  function  furnished  by  NWC 
C13  3.   This  transfer  function,  which  was  deduced  from  the 
test  data,  closely  approaches  the  actual  measurements  and  is 
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given  as 

T 
3    _  8 

THRUST     (1. 233s*l ) (50. 76s+l ) 
(In  this  expression,  thrust  in  Newtons  gives  a  temperature 
response  in  degrees  Fahrenheit  above  ambient. ) 
Vane  Thermal  Response 

Figure  4  shows  the  computed  vane  temperature  history 
using  "current  best  estimates   of  the  problem  parameters 
described  previously.   This  fi    -e  shows  that  the  predicted 
temperature-time  history  at  the  reference  shaft  node  is  well 
in  excess  of  that  obtained  in  the  tests. 

In  the  model  as  it  is  presently  configured,  the  quanti- 
ties encumbered  by  the  greatest  amount  of  uncertainty  are 
those  associated  with  the  boundaries  of  the  thermal  system 
--  the  convective  inputs  and  the  conductive  cooling  of  the 
vane  by  the  thermal  sink  effect  of  the  mount.   Some  prelimi- 
nary experimentation  with  the  simulation  has  led  to  the 
insight  that  the  chief  control  over  the  maximum  shaft 
temperature  is  the  mount  thermal  resistance.   On  the  other 
hand,  the  dynamic  temperature  rise  is  mainly  affected  by  the 
convective  thermal  resistances  acting  at  the  vane  surface. 
These  observations  lead  to  the  premise  that  the  effective 
time  constants  for  these  two  processes  are  widely  separated 
--  a  clue  to  which  is  given  in  the  previous  transfer-func- 
tion model. 

Adjustment  of  these  factors  leads  results  of  the  sort 
shown  in  Fig.  5.  which  illustrates  the  sensitivity  of  the 
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shaft  temperature  response  to  the  convective  film  coeffi- 
cients.   It  will  be  seen  from  this  figure  that  in  order  to 
obtain  reasonable  agreement  between  theory  and  experiment  it 
is  necessary  to  reduce  these  coefficients  (with  correspond- 
ing increases  in  the  respective  thermal  resistances)  by  a 

factor  on  the  order  of  60%  to  80*/.  --  the  ratio  F   in  Fig.  5. 

R 

The  need  for  this  offset  is  largely  attributable  to  omission 

of  the  cooling  effects  of  radiation  and  ablation.    It  is 

apparent,  however,  that  the  model  is  capable  of  reproducing 

the  main  transient  features  (time  constants)  of  the  vane 

thermal  response. 

Figure  6  illustrates  the  sort  of  agreement  with  test 

data  that  is  achievable  using  the  present  model.   To  achieve 

the  results  shown  in  Figs.  5  and  6, result,  a  mount  thermal 

resistance  of  R^_  =  10.3  K/W  was  used  --  a  value  not  incon- 

3G 

sistent  with  the  mass  and  configuration  of  the  mount.   Whe- 
ther or  not  the  reduction  factor  of  72%  is  consistent  with 
the  effect  of  radiative  and  ablative  cooling  is  yet  to  be 
determined.   In  any  case,  it  is  noteworthy  that  the  agree- 
ment illustrated  in  Fig.  6  is  attainable  by  means  of  a  sin- 
gle constant  factor. 

In  addition  to  providing  the  framework  of  a  working 
computational  code,  the  results  of  the  simulation  indicate 
that  it  may  be  feasible  to  predict  the  bulk  thermal  behavior 
of  selected  critical  elements  of  the  jet  vane.   In  order  to 
obtain  these  preliminary  results,  no  particular  attention 
has  been  paid  to  the  precision  with  which  the  boundary  ther- 
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mal  resistances  are  calculated.   For  the  purposes  of  this 
work,  the  simplicity  of  the  model  gives  some  encouragement 
that  it  may  be  adaptable  to  the  deductive  or  system-identi- 
fication mode  of  analysis.   Preliminary  results,  reported  in 
subsequent  sections  of  this  report,  indicate  that  this  is  in 
deed  the  case. 

VI.    STATE  SPACE  MODEL  FORMULATION 
The  results  presented  thusfar  have  been  obtained  from  a 
sequential  integration  of  the  energy  balance  relationships 
at  each  calculation  node,  and  Eq.  (2)  illustrates  the  form 
in  which  these  relationships  have  been  cast  in  the  DSL  pro- 
gram for  purposes  of  simulation.   If  the  physical  quantities 
(resistances  and  capacitances)  appearing  in  the  governing 
system  of  equations  are  truly  constant,  then  the  system  is 
linear  and  the  efficiency  of  a  state  space  formulation  may 
be  brought  to  bear.   The  matrix  formats  that  characterize 
the  state  space  formulation  are  particularly  useful  in  pro- 
viding a  systematic  presentation  of  the  mathematical  model. 
In  addition,  the  following  results  have  been  obtained 
using  the  personal  computer  version  of  the  software  package 
called  "MATRIXx"  (product  of  Integrated  Systems  Incorpor- 
ated, Palo  Alto,  Calif. ).   Incorporated  in  MATRIXx  is  a  vast 
variety  of  matrix  manipulation  schemes  that  are  often  found 
to  be  essential  in  the  modeling  and  simulation  of  dynamic 
systems  (particularly  control  systems).   Included  are  means 
for  conversion  from  continuous  to  discrete  systems  and  for 
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transformation  of  systems  fx^om  state  space  to  transfer  func- 
tion models  and  vice-versa.   The  system  identification  capa- 
bilities of  MATRIXx,  essential  to  the  goals  of  this  study, 
are  discussed  in  a  subsequent  chapter. 

Mathematical  relationships  are  expressed  in  state  space 
form  by  arranging  them  such  that  the  highest  order  deriva- 
tives of  the  dependent  variables  (the  temperatures  at  the 
various  nodes,  in  this  case)  are  given  as  linear  functions 
of  the  lower  order  derivatives.   For  node  2,  for  instance, 
Eq.   ( 1 )  may  be  written 
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C2   RF2 


)    ♦ 


R2 


2  12     2  23 


C^R_ 
2  F2 
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Referring  to  Fig.  3,  the  corresponding  relationships  for 
nodes  1  and  3  are  obtained: 
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The  internal  node  designated  (i)  in  Fig.  3,  which  was 
included  to  allow  separate  estimation  of  the  tip  and  vane 
resistances,  is  not  a  storage  node  and,  as  previously  noted, 
energy  balance  for  this  node  does  not  include  a  rate  term: 

T .  -  T- 


Rii 


(6) 


12 


Equation  (6)  may  be  solved  for  T   and  Eqs.  (3)  and  (4)  may 
be  used  to  obtain  an  expression  for  its  rate  of  change  in 
the  format  of  Eqs.   (3)-(5).   This  would  lead  to  a  four-state 
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formulation  for  the  four  node  temperatures.   For  efficiency 
of  computation,  however,  it  seems  good  practice  to  eliminate 
such  non-storage  nodes  from  consideration  --  temperatures  at 
such  nodes  can  always  be  calculated  from  algebraic  expres- 
sions, such  as  Eq.   (6),  involving  the  temperatures  at  adja- 
cent nodes.   Using  Eq.   (6)  to  eliminate  T   from  Eqs.   (3)  and 
(4)  gives  the  following  3-state  set: 

T  T         T 

4-         i#_l       1  ,      2   1       Rl   1 

1  Cl  RF1    R12      Cl  R12     Cl  RF1 

T  T  T  T 

*       11  2  1       1       1  .  3  R2  1 

T   =  -  ( ♦  ♦  )  «■  *  (8) 

■^    \-*  ,-*    R .  J-*  w  .-*  R  — -»     re.  _     R,-*^.  Lrf-_r\-.-_  Lrf  _  R__ 

2  12  2  F2     12     23  2  23  2  F2 

T           T  T 

*      2  1       3,1  1  .  G  1 

T   =  -   ( ♦  )  ♦  (9) 

3  C_  R_,_      t*— .   R  — ,»  K_  _  G_  R^« 

3   23      3    23     3G      3   3G 

where  R  _  =  R    *  Ri2* 

Arranged  in  such  a  fashion,  it  is  apparent  that  the 

nodal  expressions  possess  a  certain  amount  of  symmetry. 

Each  of  the  coefficients  of  temperature  on  the  righthand 

sides  of  Eqs.  (7)- (9)  have  the  dimension  of  inverse  time 

and,  in  fact,  the  RC  products  are  representations  of  the 

time  constants  describing  the  energy  transport  processes 

occuring  at  and  around  the  nodes.   To  further  illustrate 

this  foraal  arrangement,  it  is  useful  to  define  the 

following  parameters  (inverse  time  constants): 

a   =-(a   ♦a   +b   ).  a   =1/C  R   .  a   =0 

11     12   13   11  '   12     1  12'  13 

a21  =  1/C2R12'  a22=-<a2ra23^b22>'     a23=1/C2R23  <10) 

a31  =  °'  a32=1NC3R23'  a33=-<a3ra32  +  b33> 
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where 

bll=  1/C1RF1  '  b22  =  1/C2RF2  '  b33  =  1/C3R3G 
With  these  definitions,  Eqs.   (7)-(9)  may  be  written 

T  =  A  T  *  B  U  (11) 

where  the  "state"  consists  of  the  n  =  3  nodal  temperatures: 

T  =  (T1   T2   T3> '   , 

the  "input"  vector  contains  the  m  =  3  boundary  temperatures: 

U  =  ,TR1   TR2   V    ' 
and  the  nxn  matrix  A  and  the  nxm  matrix  B  are  given  by 


A  = 


a     a     a 
11    12    13 

a     a     a 
21    22    23 

a     a     a 
31    32    33 


B  = 


11 

0 


22 

0    b 


0 
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Equation  (11)  is  expressed  in  the  state-space  form  familiar 

to  the  analysis  of  controls  systems,  and  if  the  vector  of 

input  temperatures  U  are  specified  the  time  histories  of  the 

nodal  temperatures  are  easily  obtained. 

Because  the  "ground"  temperature  T_  is  a  specified 

u 

constant  in  this  model  (in  the  case  at  hand  it  is  290K )  the 
calculation  may  be  simplified  by  referring  all  temperatures 
to  ground.  This  has  the  effect  of  removing  the  ground  tem- 
perature as  an  input  (it  is  zero  relative  to  itself)  and  the 

input  vector  becomes  U  =  (T  .   T  _>'.   With  only  the  recov- 
f  —      rl    r2  7 

ery  temperatures  as  inputs  (m  =  2)  the  last  column  in  the 
input  coefficient  matrix  B  is  redundant  and  must  be  deleted. 

The  three  nodal  temperatures  are  "outputs, "  in  a  con- 
trol systems  sense,  and  the  usual  output  relationship  may  be 
written  as 
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y_  =  C  T  t  D  U 
where  C  is  a  3x3  identity  matrix  and  D  is  a  3x2  null  matrix. 
As  a  final  step  prior  to  simulation,  the  four  matrices  may 
be  combined  into  a  single  partitioned  matrix,  S,  containing 
the  four  submatrices  A,  B_,  C,  and  D.   This  is  particularly 
useful  when  dealing  with  the  MATRIXx  software.   Thus 


S  = 


A  '  B 


C  ,  D 


or 


S  = 


a     a     a     b 
11    12    13    11 


21 

a22 

B23 

0 

b22 

31 

a32 

a33 

0 

0 

1 

0 

0 

0 

0 

0 

1 

0 

0 

0 

0 

0 

1 

0 

0 

(12) 


In  this  formulation  the  thermal  characteristics  of  the 

system  are  specified  in  the  A  matrix  and  the  input  resis- 
tances are  found  in  the  elements  making  up  the  B  matrix. 
Note  that  the  diagonal  elements  of  the  A  matrix  are  the 
negative  inverses  of  the  characteristic  time  constants  for 
the  thermal  response  of  each  of  the  nodes.   As  defined  in 
Eq.  (10),  these  diagonal  elements  are  sums  of  the  values  on 
the  corresponding  row,  including  the  input  matrix.   This 
relationship  expresses  the  conservation  of  energy  embodied 
in  the  analysis  and  the  diagonal  elements  cannot  be  indepen- 
dently varied  without  violating  this  principle.   This  is  an 
important  constraint  when  system  identification  procedures 
are  applied  to  the  system. 


34 


The  generality  of  the  formulation  is  apparent,  and 
there  would  be  no  particular  difficulty  in  adding  storage 
nodes  to  the  model  in  order  to  obtain  a  more-detailed  view 
of  the  temperature  distribution  within  the  vane.   Such 
refinement  at  this  point  seems  unwarranted,  however,  since 
little  information  is  available  to  allow  a  concomittant 
refinement  of  the  distribution  of  the  thermal  inputs.    In 
fact,  the  present  model  would  not  be  drastically  affected  if 
the  stagnation  point  and  boundary  layer  heat  transfer  were 
combined  into  a  single  input,  thus  reducing  the  model  to  a 
two-state  system. 


SIMULATION  IN  STATE  SPACE 
Figure  3  includes  a  table  of  the  thermal  resistances 
and  nodal  capacitances  that  have  been  estimated  for  the 
problem  at  hand.    (In  estimating  the  "effective"  input 
resistances  given  in  Fig.  3,  the  calculated  heat  transfer 
coefficients  have  been  reduced  by  the  72%  factor  leading  to 
the  results  described  in  the  previous  chapter. )   When  these 
values  are  incorporated  into  the  present  analysis,  the  A  and 
B  matrices  take  the  following  values: 


A  = 


-5.15    4.7      0 
0. 23   -0. 61    0. 18 
0      0.38   -0.48 


B  = 


0.  45    0 
0     0.  20 
0      0 


(  13) 


For  the  tip  node  the  input  temperature  is  the  stagna- 
tion temperature  ( 2650K )  and  for  the  vane  the  recovery  tem- 
perature calculated  by  the  methods  of  Ch.  IV  is  2550K  (very 
nearly  the  same  as  the  stagnation  temperature  in  this  high- 
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speed  compressible  flow).  When  referred  to  an  ambient  tem- 
perature of  290K,  the  appropriate  input  values  for  tempera- 
ture are 

T_,  =  2360K    and   T__  =  2260K 
Although  the  thrust  was  used  as  the  fundamental  input  for 
the  simulations  previously  described,  these  temperatures  may 
be  assumed  to  follow  an  input  schedule  that  is  more-or-less 
synchronous  with  the  thrust.   Accordingly,  in  the  simula- 
tions that  follow  the  recovery  temperatures  have  been  given 
the  same  ramp-up,  plateau  (at  the  values  given  above),  and 
ramp-down  contour  that  is  described  in  Ch.  V. 
Linear  Time-Invariant  (LTI)  Model 

With  the  model  matrices  values  given  in  Eq.  (13)  and 
the  scheduled  input  temperatures  described  above,  simulation 
of  the  vane  thermal  response  is  accomplished  using  the  MAT- 
RIXx  LSIM  command  C 14 ,  p. 10-63.   In  this  case,  the  syntax  is 

[T,Y]  =  LSIM(S,NS,  U,  DELTAT) 
where  S  is  the  system  matrix  formed  as  in  Eq.  (12)  and  in- 
corporating the  values  given  in  Eq.  (13);  NS  is  the  number 
of  states  ( NS  =  3  --  the  three  nodal  temperatures);  U  is  the 
input  array,  in  this  case  expressing  the  scheduled  input 
temperatures;  and  DELTAT  is  the  specified  time  increment 
between  steps  in  the  U  array.   In  all  simulation  cases  de- 
scribed herein,  a  time  increment  of  0. Is  has  been  used  for  a 
duration  of  15s.   Thus  the  time  vector  has  dimensions  151x1 
and  the  U  array  is  151x2  for  the  two  inputs. 

Figure  7  shows  the  results  of  the  simulation  assuming 
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that  the  elements  of  the  A  and  B  matrices  are  constant  with 
time.   As  expected,  the  results  are  unsatisfactory  because 
the  linear  time-invariant  (LTI)  assumption  does  not  account 
for  the  variation  of  the  heat  transfer  coefficients  with 
time.    In  other  words,  the  elements  of  the  B  matrix  and  the 
diagonal  elements  of  the  A  matrix  must  be  calculated  at  each 
time  step  in  the  simulation. 
Linear  Time-Variant  (LTV)  Model 

With  MATRIXx,  the  accomodation  of  time-varying  parame- 
ters may  be  relatively  easily  achieved  by  using  the  SYSTEM 
BUILD  capability.   SYSTEM  BUILD  provides  an  interactive, 
menu-driven  graphical  environment  for  building,  modifying, 
and  editing  computer  simulation  models.   Any  combination  of 
linear,  non-linear,  continuous-time,  discrete-time,  or  mul- 
ti-rate models  that  describe  a  system  may  be  constructed 
from  a  library  of  basic  building  blocks.   Systems  may  be 
modeled  by  dividing  them  into  individual  components  with 
each  component  being  described  by  a  specific  type  of  func- 
tional block.   Super-blocks  can  be  used  to  represent  as- 
semblies of  individual  blocks  and  the  heirarchical  design 
provides  an  organized  representation  of  the  physical  system 
that  is  modeled. 

In  SYSTEM  BUILD,  any  super-block  or  set  of  nested  su- 
per-blocks can  be  simulated,  linearized,  and  analyzed.   No 
user  coding  is  necessary  (although  it  is  allowed),  and  SYS- 
TEM BUILD  can  be  driven  by  user-defined  command  procedures 
that  automatically  build  or  modify  system  models.   All  of 
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this  capability  is  available  for  PC  hardware  (which  has  been 
used  thusfar  in  this  study)  and  the  net  effect  of  the  SYSTEM 
BUILD  capability  is  that  the  user  may  concentrate  efforts  on 
design,  analysis,  and  simulation  in  an  efficient  graphical 
environment . 

The  SYSTEM  BUILD  Model.    In  this  study,  the  vane  ther- 
mal behavior  has  been  modeled  by  means  of  seven  super-blocks 
named  as  follows: 

N0D1IN,  N0D2IN,  N0D3IN 
N0DE1,  N0DE2,  N0DE3 
VANE5 
The  appendix  provides  illustrations  of  these  super-blocks, 
and  a  brief  discussion  will  be  provided  here. 

The  first  three  blocks  are  compute  the  parameters  ne- 
cessary for  the  state  equation.   The  input  to  both  N0D1IN 
and  N0D2IN  is  the  ramp-up,  plateau,  ramp-down  schedule  pre- 
viously described,  but  with  a  plateau  value  of  unity.   In 
N0D1IN,  for  example,  this  input  is  multiplied  by  T    and  b 

HI         11 

as  gains  to  form  the  time-varying  element  in  the  first  row 
of  the  input  coefficient  array  B.   Also  in  N0D1IN,  the  quan- 
tity a  _  is  generated  by  means  of  the  step  function  and 
this,  in  turn,  is  used  to  form  -a    =  a  _♦  b    as  given  in 
Eq.   (10).   The  outputs  of  N0D1IN  are  -a..  ,  b   T   ,  and  a   . 

11      11  Kl  L  jl 

N0D2IN  and  N0D3IN  are  similar  in  function,  with  the  excep- 
tion that  N0D3IN  does  not  require  the  time-varying  input. 

The  first  three  super-blocks  listed  above  are  elements 
of  the  second  three,  respectively.   In  N0DE1,  for  example, 
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the  external  input  T   is  combined  with  the  outputs  of  NQD1IN 

to  form  the  quantity  a^7^  ♦  a12T2  *  bllTRl  =  "*V   An  inte" 

gration  is  also  within  N0DE1  to  obtain  T,  from  T.  .   T.  is 

111 

the  output  of  the  N0DE1  super-block  and,  in  a  similar  way, 
NQDE2  and  N0DE3  provide  outputs  T   and  T  . 

Super-block  VANE5  "wires"  N0DE1,  N0DE2,  and  N0DE3  to- 
gether in  order  to  provide  simultaneous  solution  for  the 
three  nodal  temperatures.    In  VANE5  (see  appendix),  NQDE1 
receives  the  dual  ramp  input  and  T   while  generating  T  ; 
N0DE2  receives  the  dual  ramp  input,  T  ,  and  T   while  gener- 
ating T  ;  and  N0DE3  receives  T   and  outputs  T  . 

To  execute  the  SYSTEM  BUILD  model,  the  program  is  first 
checked  by  means  of  the  ANALYZE  command  C l5 ] .   Following 
this,  one  of  six  integration  methods  may  be  selected  from 
the  menu  provided.   The  results  described  here  have  been 
obtained  with  a  variable-step  Kutta-Merson  scheme  which  is 
an  explicit  fourth-order  method  that  employs  the  largest 
time  step  possible  while  remaining  within  error  tolerances 
(maximum  step  size  is  the  time  increment  specified  for  the 
simulation ) . 

Simulation  is  achieved  by  executing  the  SYSTEM  BUILD 
command  SIM,  as  in  "Y  =  SIM(T,U)",  where  the  time  vector  T 
and  the  input  U  must  be  provided.   In  the  case  at  hand,  T  is 
the  151x1  vector  containing  15  seconds  of  0.1-second  steps 
and  U  is  the  151x2  array  exhibiting  the  ramp-up,  unity  pla- 
teau, ramp  down  schedule  previously  described. 

Figure  8  shows  the  results  of  simulation  using  the 
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SYSTEM  BUILD  model  to  provide  linear  time-variant  (LTV) 
behavior  of  the  vane  heat  transfer  coefficients.    In  this 
figure  the  improvements  attained  with  the  SYSTEM  BUILD  LTV 
model  are  apparent.   With  the  model  now  cast  in  this  format 
it  is  possible  to  take  the  system  identification  approach 
that  has  been  the  ultimate  goal  of  this  study. 

VIII.    PARAMETRIC  SYSTEM  IDENTIFICATION 
The  philosophy  of  parameteric  system  identification 
(PSI)  is  to  use  measured  response  data  to  "massage"  the 
parameters  of  a  system  model  such  that  simulation  results 
match  those  of  experiments.   Some  of  this  has  already  been 
accomplished  by  the  manual  variation  of  parameters  described 
in  Ch.  V,  but  PSI  will  now  be  used  to  complete  this  process 
in  a  systematic  way. 

In  the  problem  at  hand,  it  is  reckoned  that  the  most 
uncertain  parameter  values  are  those  associated  with  the 
heat  transfer  processes  at  the  vane  surface  and  the  thermal 
sink  effects  of  the  vane  mount.   The  vane  temperature  to  be 
matched  is  that  of  the  shaft  <T  ),  and  although  the  stagna- 
tion heat,  transfer  process  at  the  vane  tip  has  been  modeled 
this  provides  relatively  little  heat  flux  into  the  vane  and 
is  therefore  of  secondary  importance  as  far  as  the  thermal 
response  at  the  shaft  is  concerned. 

Thus  the  results  described  here  have  been  obtained  by  a 
two-parameter  variation  involving  b    and  a   .   It  is  to  be 

recalled  that  the  variation  of  b__  will  also  require  a 

22  M 
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change  in  a~.-,  in  accordance  with  the  energy  balance,  and 
that  this  is  accomodated  in  the  SYSTEM  BUILD  model.    In  ad- 
dition, the  value  of  a    has  within  it  b    which  is  not  used 
in  the  simulation  since  the  relative  ground  temperature  is 
zero.   Nevertheless,  a  knowledge  of  a    will  lead  to  the 


appropriate  value  of  b    since  b 


a33  *  Q32  aS  9±Ven  ln 


Eq.   (10). 

The  MATRIXx/PC  family  includes  a  variety  of  capabili- 
ties that  provide  for  data  analysis,  parametric  and  non-par- 
metric  system  identification,  and  digital  filtering.   In 
this  work,  PSI  is  achieved  by  using  a  scheme  known  as  "maxi- 
mum likelihood  estimation"  (MLE)  Cl6],   As  its  name  implies, 
MLE  maximizes  the  likelihood  of  the  parameter  estimates 
given  the  observations  (in  this  case,  the  measured  transient 
response  of  T  ).   MLE  operates  on  an  entire  time  history 
(rather  than  one  sample  at  a  time)  and  is  applicable  to 
multiple-input /multiple-output  (MIMQ)  systems  such  as  that 
under  analysis  here. 

The  command  that  executes  the  MATRIXx  version  of  the 
MLE  procedure  is  MAXLIKE  [17],  with  the  following  syntax: 

[YHAT,P]  =  MAXLIKE(U, T3, P0, ' VBM2', NIT) 
Items  to  be  provided  by  the  user  are  U,  containing  the  time 
vector  and  array  of  input  coefficients;  T3,  the  vector  of 
observed  values  that  are  the  standards  for  comparison;  P0, 
the  vector  of  initial  estimates  for  the  parameters  to  be 
varied ;' VBM2 ' ,  a  user-supplied  command  file  that  creates  the 
values  of  T3  based  upon  U  and  the  current  values  of  P;  and 
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NIT,  the  specified  number  of  iterations  to  be  performed. 

The  results  of  MAXLIKE  are  YHAT,  the  system  response 
using  current  values  of  the  parameters;  P,  which  contains 
these  current  values;  and  RSS  which  is  the  sum  of  the 

squares  of  the  differences  between  the  measured  and  pre- 

2 
dieted  values  at  each  time  step  --  RSS  =  sum(T3  -  YHAT)  . 

As  has  been  shown,  the  LTV  nature  of  the  vane  heat 
transfer  problem  requires  the  use  of  SYSTEM  BUILD.   Thus,  in 
the  problem  at  hand,  the  command  file  'VBM2'  reads  as  fol- 
lows : 

BUILD, EDIT, N0D2IN, 0, EXAMINE, 2, NEW,  P( 1 ),  DONE,  TOP,  .  .  . 
EDIT, N0D3IN, 0, EXAMINE, 3, NEW, P(2), DONE,  TOP,  .  .  . 
ANALYZE, VANES 

Y3=SIM(U( :, 1 ), U< :, 2:3) ) ; 

YHAT=Y3( :, 3) ; 

RETURN 
Although  this  language  may  not  be  familiar  to  the  reader,  it 
may  be  possible  to  deduce  that  in  the  SYSTEM  BUILD  lexicon 
these  instructions  sequentially  update  P<1)  (=  b   )  which  is 
in  position  2  of  super-block  N0D2IN  and  P(2)  (=  a~3>  which 
is  in  position  3  of  super-block  N0D3IN.   Y3  is  the  name  gi- 
ven to  the  three-temperature  array  obtained  from  simulation 
using  the  current  parameter  values,  P(l)  and  p(2),  with  YHAT 
(=  T  )  contained  in  third  column.   It  should  be  noted  that 
in  each  iteration  the  MAXLIKE  procedure  executes  'VBM2'  at 
once  for  each  parameter  value  under  consideration. 

Figure  9  is  the  result  of  a  PSI  run  of  three  iterations 
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using  initial  parameter  values  given  previously  ( b    =  0.20 
and  a  _  =  -0.48).   The  "maximum  likelihood"  values  obtained 
through  the  use  of  PSI,  and  giving  the  results  shown  in  Fig. 
9,  are  b    =  0. 2027  and  a    =  -0. 4565.   The  value  of  the 
residual  associated  with  these  results  is  RSS  =  59.4.    (This 
is  a  very  small  value,  bearing  in  mind  that  it  is  a  sum  of 
151  squared  temperature  differences. ) 
Discussion  of  the  PSI  Results 

As  was  done  in  Ch.  V,  the  experimental  shaft  tempera- 
ture data  (shown  as  x's  in  Figs.  7-9)  have  been  obtained 
from  the  transfer  function  provided  by  NWC.   In  this  case, 
the  temperature  above  ambient,  in  Kelvins,  is  given  by 

T_ 


4.  48 


(14) 


THRUST    (1.  233s+l)  (50.  76s«-l) 
where,  as  before,  the  thrust  has  the  dimension  of  Newtons 
(see  Fig.  6).   A  state-space  equivalent  of  this  expression, 
which  is  easily  obtained  by  means  of  the  MATRIXx  function 
SF0RM,  is 

x  =  Ax  +  Bu 

y  =  Cx  ♦  Du 
where 


A  = 


-0.8307  -.0.0160 
1  0 


B  = 


1 
0 


C  =  CO    0.0716] 


and  D  =  0.   (In  this  representation,  y  is  the  shaft  tempera- 
ture T   and  u  is  the  thrust. )   In  the  combined  matrix  format 
previously  discussed  in  connection  with  Eq.  (12),  these  ma- 
trices yield 
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S  = 


0.8307  -0.0160  1 
1  0  0 
0     0.0716   0 


(  15) 


The  agreement  illustrated  in  Fig.  9,  which  is  clearly 
excellent,  can  be  evaluated  in  several  ways.   For  instance, 
the  temperature  response  resulting  from  the  identified  sys- 
tem may  be  cast  in  a  form  that  is  directly  comparable  with 
that  provided  by  NWC.   To  accomplish  this,  the  system  model 
may  be  converted  to  the  2-state  matrix  form  given  above,  and 
the  three  parameters  of  S  subjected  to  system  identifica- 
tion.  This  has  been  done  by  means  of  the  MAXLIKE  procedure 
with  the  following  result: 


S'  = 


-0.8218  -0.0173  1 
1  0  0 
0     0.0716   0 


with  a  value  of  RSS  =  63. 4.   Comparison  of  this  result  with 
Eq.   (15)  demonstrates  the  agreement.   Using  the  TFORM  com- 
mand of  MATRIXx,  the  state-space  representation  may  be  con- 
verted back  to  transfer-function  form,  yielding 


V 


4.  13 


THRUST    ( 1. 250s+l ) (46. 08s+l) 
Comparing  with  Eq.  (14)  the  agreement  is,  again,  remarkable. 
The  reader  should  bear  in  mind  that  these  arguments  are  made 
here  only  for  the  purpose  of  demonstrating  the  effectiveness 
of  the  PSI  procedure  in  this  case.   The  various  transfer 
functions  and  parameters  in  these  models  have  no  relation  to 
the  physical  situation  --  in  fact,  since  they  represent  LTI 
syetems,  they  are  wholly  inadequate  for  any  purpose  other 
than  expressing  the  experimental  measurements  simulation 
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results  in  functional  forms. 

Several  sets  of  initial  values  have  been  used  in  iden- 
tifying the  parameter  values  in  the  physical  model,  and  the 
results  have  been  found  to  be  independent  of  this  factor. 
In  addition,  a  number  of  tests  have  been  performed  to  deter- 
mine the  effect  of  varying  other  parameters  in  the  system 
model  and  no  significant  improvement  has  been  obtained  on 
the  results  illustrated  in  Fig.  9.   The  selection  of  the 
appropriate  parameters  to  vary  in  the  PSI  process  requires 
some  appreciation  for  the  physical  process  modeled.   Thus  it 
has  been  observed,  for  instance,  that  the  selection  of  too 
many  parameters  will  lead  to  unrealistic  values,  such  as 
negative  thermal  conductances,  unless  constraints  are  added 
to  the  MAXLIKE  procedure.   These  and  other  aspects  of  the 
method  are  under  continued  study. 

IX.    CONCLUSIONS  AND  RECOMMENDATIONS 
The  emphasis  in  the  work  reported  here  has  been  to  in- 
vestigate the  feasibility  of  using  system  identification 
methods  to  develop  tools  for  use  in  the  design  of  jet  vane 
TVC  devices.   Before  proceeding  to  conclusions  in  this 
regard,  however,  it  may  be  useful  to  comment  on  the  results 
obtained  in  the  case  that  has  been  analyzed. 
Results  of  the  Analysis 

The  excellent  agreement  between  measured  and  predicted 
shaft  temperature  responses  shown  in  Fig.  9  gives  some  cre- 
dibility to  the  model  parameters  used  in  the  simulation. 
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These  may  in  turn  be  related  to  vane  physical  properties  to 
lend  some  insights  regarding  those  quantities  that  were 
initially  subject  to  relatively  large  uncertainties.   For 
the  heat  sink  effect  of  the  shaft  mount,  for  instance,  the 
parameter  b    =  -  a    -  a~-  may  be  obtained  from  the  identi- 
fied value  of    a    =  -0.46  and  the  assumed  (and  relatively 
certain)  value  of  a    =  0.38.    The  result  is  b    ~    0.08  s 
and  this  may  in  turn  be  related  to  the  vane  physical  proper- 
ties since  b„_  =  1/C_R__.   This  new  value  is  down  from  the 
J3         3  Ju 

original  value  of  0.10  s    (Fig.  3),  and,  since  RC  =  mcL/kA, 
the  result  may  be  used  to  adjust  the  rationale  used  in  esti- 
mating the  effective  mass  and  geometry  of  the  thermal  sink 
representing  the  mount  (assuming  that  the  heat  capacity  and 
thermal  conductivity  are  firm  values). 

The  other  parameter  identified  in  the  PSI  analysis  is 
the  vane  thermal  input  coefficient,  b    =  0. 2027  s   .    If 

this  value  is  used  to  recalculate  the  convective  heat  trans- 

2 

fer  coefficient,  a  value  of  h  =  638  W/m  K  is  obtained  which 

is  about  72%  less  than  that  estimated  in  Ch.  IV  for  turbu- 
lent boundary  layer  convection.   The  conclusion,  as  has  been 
mentioned,  is  that  the  complex  processes  of  ablation  and 
radiation  account  for  a  considerable  amount  of  cooling  of 
the  vane.    It  is  important  to  keep  in  mind  that  the  single 
values  quoted  here  are  useful  only  inasmuch  as  they  give 
good  results  for  predicting  shaft  temperature  response  using 
the  model  that  has  been  constructed.    In  this  case,  the  mo- 
del includes  a  ramp-up,  plateau,  ramp-down  behavior  for  b 
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so  that  the  above  result,  which  applies  to  the  plateau  value 
of  the  "effective"  heat  transfer  coefficient,  accounts  in  a 
crude  way  for  the  fact  that  vane  heat  transfer  effects  vary 
widely  throughout  the  firing.   The  assumption  of  a  constant 
plateau  value  for  h  is,  of  course,  subject  to  question  and 
refinement  if  sufficient  information  is  available. 

Referring  to  the  vane  surface  temperature  responses 
shown  in  Fig.  9,  it  is  seen  that  a  maximum  of  about  1240K 
( 950K  above  ambient)  is  reached  at  about  3.5  seconds  (3  se- 
conds after  motor  ignition).   This  value  is  well  below  the 
melting  point  of  the  vane  material,  and  the  tensile  modulus 
is  still  about  22,000  ksi  at  this  temperature  (though  it  is 
well  below  the  ambient  value  of  about  45,000  ksi).    It  might 
be  predicted  from  these  results  that  the  1/4-scale  vane  that 
has  been  modeled  would  remain  more-or-less  intact  during  a 
firing  of  the  type  that  has  been  simulated.   A  noteworthy 
feature  of  these  developments  is  that  it  would  be  a  rela- 
tively simple  matter  to  simulate  the  thermal  response  of  a 
full-scale  vane.   Each  of  the  parameters  in  the  present  mo- 
del is  amenable  to  direct  scaling  procedures  such  as  those 
that  have  been  given  in  Ch.  IV. 
The  PSI  Method 

The  results  thusfar  obtained  in  this  study  justify  some 
confidence  that  the  thermal  behavior  of  the  jet  vane  can  be 
modeled  using  a  relatively  simple  and  straight-forward  model 
structure.   But  although  the  structure  of  the  model  is  well 
defined,  some  of  the  values  of  the  system  parameters  neces- 
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sary  to  make  the  model  "work"  are  subject  to  considerable 
uncertainty.  In  this  study  the  power  of  the  system  identi- 
fication method  has  been  brought  to  bear  in  removing  some  of 
this  uncertainty.  In  addition,  it  is  worth  mentioning  that 
the  relatively  new  capabilities  now  available  with  software 
products  like  MATRIXx  have  greatly  facilitated  these  advan- 
ces. 

Future  research  in  these  areas  is  warranted  from  both 
the  modeling  and  system  identification  points  of  view.   Us- 
ing the  model  that  has  been  developed  here,  or  perhaps  an 
even  simpler  model  (without  the  tip  node),  it  is  important 
to  determine  if  the  good  results  obtained  in  the  present 
instance  can  be  expected  in  other  cases.   With  this  in  mind, 
the  results  of  other  firing  tests  should  be  predicted  with 
the  model  and  further  system  identification  conducted  if  it 
proves  to  be  necessary.   The  basic  issue  would  be  to  deter- 
mine the  robustness  (or  lack  thereof)  of  the  present  model. 

From  the  modeling  point  of  view,  a  number  of  inter- 
esting questions  remain  of  a  more  theoretical  nature.   The 
relationship  between  the  complexity  of  the  model  and  the 
adequacy  of  experimental  comparison  data  is  of  some  impor- 
tance.  In  addition,  more  realism  could  be  introduced  into 
the  present  model  in  an  effort  to  determine  the  sensitivity 
of  the  results  to  such  matters.   For  instance,  it  might  be 
possible  to  introduce  a  multiple  thermal  input  model  that 
accounts  for  radiation  and  ablation  as  well  as  thermal  con- 
vection.  Other  elements  of  realism  could  be  introduced  by 
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allowing  for  the  temperature  dependency  of  several  of  the 
parameters  now  taken  as  constant.   These  include  the  thermal 
conductivity  of  the  vane  material  and  the  reference  tempera- 
ture used  in  calculating  the  thermal  properties  of  the  motor 
exhaust  gases. 

Another  aspect  of  the  method  that  may  lead  to  important 
insights  is  that  the  initial  and  final  thermal  responses 
seem  to  he  somewhat  uncoupled  under  the  present  circumstan- 
ces.  For  example,  the  LTI  response  illustrated  in  Fig.  7  is 
quite  adequate  during  the  initial  phase.   The  implication  is 
that  the  initial  and  final  transients  might  be  used  to  iden- 
tify particular  parameter  values  under  the  assumption  that 
they  are  time-invariant.    If  this  were  true,  the  handling  of 
nonlinear  behaviors,  should  this  become  necessary,  would  be 
greatly  facilitated. 
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APPENDIX 
Diagrams  of  System  Build  Super-Blocks 
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